%     VIRTUAL FABRIC - программа для получения эффективных свойств текстильных композитов
%     Copyright (C) 2020  Михаил Козлов
%     Сайт проекта: <http://www.fastenerdoctor.ru/>
%     Адрес электронной почты для связи с разработчиком:
%     <doctor.fastener@gmail.com>
% 
%     Эта программа является свободным программным обеспечением: Вы можете
%     распространять ее и/или изменять, соблюдая условия Универсальной общественной
%     лицензии GNU, опубликованной Фондом свободного программного обеспечения,
%     либо редакции 3, либо (на Ваше усмотрение) любой редакции, выпущенной позже.
% 
%     Эта программа распространяется в надежде на то, что она окажется полезной,
%     но БЕЗ КАКИХ-ЛИБО ГАРАНТИЙ, включая даже подразумеваемую гарантию ТОВАРНОЙ ПРИГОДНОСТИ,
%     либо ПРИГОДНОСТИ ДЛЯ ОПРЕДЕЛЕННЫХ ЦЕЛЕЙ. Ознакомьтесь с Универсальной общественной
%     лицензией GNU для получения более подробной информации.
% 
%     Вы должны были получить копию Универсальной общественной лицензии GNU
%     вместе с этой программой. Если Вы ее не получили, то перейдите по адресу:
%     <http://www.gnu.org/licenses/>.
% 
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
% 
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
%     GNU General Public License for more details.
% 
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see <https://www.gnu.org/licenses/>.

function Yarn_Properties = Calculate_Yarn_Properties(V_resin, Analysis_Type, Material_Properties, Work_Directory, Chemical_Shrinkage)

global ErrCode;

Yarn_Properties.E1 = [];
Yarn_Properties.E2 = [];
Yarn_Properties.E3 = [];

Yarn_Properties.G12 = [];
Yarn_Properties.G13 = [];
Yarn_Properties.G23 = [];

Yarn_Properties.Nu12 = [];
Yarn_Properties.Nu13 = [];
Yarn_Properties.Nu23 = [];

Yarn_Properties.alpha11 = [];
Yarn_Properties.alpha22 = [];
Yarn_Properties.alpha33 = [];

Yarn_Properties.beta11 = [];
Yarn_Properties.beta22 = [];
Yarn_Properties.beta33 = [];

Yarn_Properties.Density = [];

Yarn_Properties.Cp = [];

Yarn_Properties.Lambda_11 = [];
Yarn_Properties.Lambda_22 = [];
Yarn_Properties.Lambda_33 = [];

% ----------------------------------- Параметры ячейки периодичности ------------------------------------

Fiber_Diameter             = 5e-003;      % Диаметр волокна

Fiber_Packing              = 'HEXAGONAL'; % HEXAGONAL | SQUARE    % Модель размещения волокон (см. руководство пользователя)
RVE_Thickness_Model        = 'AUTO';      % AUTO | SPECIFY        % Толщина ячейки периодичности

    % AUTO - толщина ячейки периодичности вычисляется автоматически
    % SPECIFY - толщина ячейки периодичности явно задается пользователем

    % Если RVE_Thickness_Model == SPECIFY
    RVE_Thickness_Value     = 1e-003;    % Высота ячейки периодичности

% ------------------------------------------- Параметры сетки -------------------------------------------

Element_Order              = 'LINEAR'; % LINEAR | QUADRATIC    % Порядок конечных элементов
Mesh_Seed_Parameter        = 10; % Параметр густоты сетки (число конечных элементов по контуру целого волокна будет составлять 12*Mesh_Seed_Parameter в случае Fiber_Packing = HEXAGONAL и 8*Mesh_Seed_Parameter в случае Fiber_Packing = SQUARE)
N_Elements_High            = 8; % Рекомендуемое пользователем число конечных элементов по толщине ячейки периодичности (используется, если RVE_Thickness_Model == SPECIFY)
Element_Critical_Ratio     = 20; % Максимально допустимое соотношение сторон конечного элемента

% -------------------------------- Параметры массива ячеек периодичности --------------------------------

RVE_Array_1                = 1; % Число ячеек по направлению 1
RVE_Array_2                = 1; % Число ячеек по направлению 2
RVE_Array_3                = 1; % Число ячеек по направлению 3

% ----------------------------------------- Параметры расчета -------------------------------------------
                           
Mechanical_BC_TYPE         = 'SOFT'; % SOFT | HARD    % Тип граничных условий в механическом расчете

                           % SOFT - периодические граничные условия со свободной границей ячейки
                           % HARD - периодические граничные условия с фиксированной границей ячейки (использовать с большой осторожностью!)
                           
Constraint_Directions      = 'XYZ';    % Направления, задающие грани с фиксированной границей в случае Mechanical_BC_TYPE == HARD. Строка, состоящая из символов X, Y, Z.

% ----------------------------------------- Параметры запуска -------------------------------------------

CPU                        = 1;                   % Число используемых ядер для расчета в Abaqus

% -------------------------------------------------------------------------------------------------------

% Проверка параметров
if ((V_resin <= 0) || (V_resin >= 1))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра V_resin !\n');
    return;
end

if (Fiber_Diameter <= 0)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Fiber_Diameter !\n');
    return;
end

switch Fiber_Packing
    case {'HEXAGONAL', 'SQUARE'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Fiber_Packing !\n');
        fprintf (1, 'Разрешенные значения: HEXAGONAL, SQUARE\n');
        return;
end

switch RVE_Thickness_Model
    case {'AUTO', 'SPECIFY'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Thickness_Model !\n');
        fprintf (1, 'Разрешенные значения: AUTO, SPECIFY\n');
        return;
end

if (RVE_Thickness_Value <= 0)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Thickness_Value !\n');
    return;
end

switch Element_Order
    case {'LINEAR', 'QUADRATIC'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Element_Order !\n');
        fprintf (1, 'Разрешенные значения: LINEAR, QUADRATIC\n');
        return;
end

if (Mesh_Seed_Parameter ~= floor(Mesh_Seed_Parameter))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Mesh_Seed_Parameter !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    return;
end

if (Mesh_Seed_Parameter < 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Mesh_Seed_Parameter !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    return;
end

if (N_Elements_High ~= floor(N_Elements_High))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра N_Elements_High !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    return;
end

if (N_Elements_High < 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра N_Elements_High !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    return;
end

if (Element_Critical_Ratio <= 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Element_Critical_Ratio !\n');
    return;
end

if (RVE_Array_1 ~= floor(RVE_Array_1))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Array_1 !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    return;
end

if (RVE_Array_1 < 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Array_1 !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    return;
end

if (RVE_Array_2 ~= floor(RVE_Array_2))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Array_2 !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    return;
end

if (RVE_Array_2 < 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Array_2 !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    return;
end

if (RVE_Array_3 ~= floor(RVE_Array_3))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Array_3 !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    return;
end

if (RVE_Array_3 < 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра RVE_Array_3 !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    return;
end

switch Analysis_Type
    case {'MECHANICAL', 'HEAT_TRANSFER', 'MESH_GENERATOR'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Analysis_Type !\n');
        fprintf (1, 'Разрешенные значения: MECHANICAL, HEAT_TRANSFER, MESH_GENERATOR\n');
        return;
end

switch Mechanical_BC_TYPE
    case {'SOFT', 'HARD'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Mechanical_BC_TYPE !\n');
        fprintf (1, 'Разрешенные значения: SOFT, HARD\n');
        return;
end

if ((numel(Constraint_Directions) < 1) || (numel(Constraint_Directions) > 3))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Constraint_Directions !\n');
    fprintf (1, 'Разрешенные значения: строка, состоящая из символов X, Y, Z.\n');
    return;
end

for i=1:numel(Constraint_Directions)
    switch Constraint_Directions(i)
        case {'X', 'Y', 'Z'}
            % Допустимые значения
        otherwise
            ErrCode = -1;
            fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Constraint_Directions !\n');
            fprintf (1, 'Разрешенные значения: строка, состоящая из символов X, Y, Z.\n');
            return;
    end
end

switch Chemical_Shrinkage
    case {'TRUE', 'FALSE'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра Chemical_Shrinkage !\n');
        fprintf (1, 'Разрешенные значения: TRUE, FALSE\n');
        return;
end

if (CPU ~= floor(CPU))
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра CPU !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    return;
end

if (CPU < 1)
    ErrCode = -1;
    fprintf (1, '\nОшибка функции Calculate_Yarn_Properties: недопустимое значение параметра CPU !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    return;
end

% Создание базовых точек
Points = Construct_Geometry(V_resin, Fiber_Diameter, Fiber_Packing);

if (ErrCode < 0)
    return;
end

% Построение сетки
[NODES, ELEMENTS, RVE_Coordinates, Node_Merge_Tolerance] = Generate_Mesh(Fiber_Packing, Points, RVE_Thickness_Model, RVE_Thickness_Value, Element_Order, Mesh_Seed_Parameter, N_Elements_High, Element_Critical_Ratio, RVE_Array_1, RVE_Array_2, RVE_Array_3, Analysis_Type);

if (ErrCode < 0)
    return;
end

% Корректировка Work_Directory
if (Work_Directory(length(Work_Directory)) == '\')
    Work_Directory = Work_Directory(1:length(Work_Directory)-1);
end

% Проверка свойств материала
Check_Material_Properties(Material_Properties);

if (ErrCode < 0)
    return;
end

% Построение систем координат
Coordinate_Systems = Define_Coord_Systems();

% Ветвление по типу анализа
switch Analysis_Type
    
    case {'MECHANICAL', 'MESH_GENERATOR'}
        
        [Density, Properties_Alpha, Properties_Beta] = Mechanical_Procedure_Yarn(Work_Directory, CPU, NODES, ELEMENTS, Material_Properties, Coordinate_Systems, RVE_Coordinates, Node_Merge_Tolerance, Chemical_Shrinkage, Mechanical_BC_TYPE, Constraint_Directions, V_resin);
        
        if (ErrCode < 0)
            return;
        end
        
        Yarn_Properties.E1 = Properties_Alpha.E1;
        Yarn_Properties.E2 = Properties_Alpha.E2;
        Yarn_Properties.E3 = Properties_Alpha.E3;
        
        Yarn_Properties.G12 = Properties_Alpha.G12;
        Yarn_Properties.G13 = Properties_Alpha.G13;
        Yarn_Properties.G23 = Properties_Alpha.G23;
        
        Yarn_Properties.Nu12 = Properties_Alpha.Nu12;
        Yarn_Properties.Nu13 = Properties_Alpha.Nu13;
        Yarn_Properties.Nu23 = Properties_Alpha.Nu23;
        
        Yarn_Properties.alpha11 = Properties_Alpha.alpha11;
        Yarn_Properties.alpha22 = Properties_Alpha.alpha22;
        Yarn_Properties.alpha33 = Properties_Alpha.alpha33;
        
        if (strcmp(Chemical_Shrinkage, 'TRUE'))
            Yarn_Properties.beta11 = Properties_Beta.alpha11;
            Yarn_Properties.beta22 = Properties_Beta.alpha22;
            Yarn_Properties.beta33 = Properties_Beta.alpha33;
        end
        
        Yarn_Properties.Density = Density;
        
    case 'HEAT_TRANSFER'
        
        [Density, Cp, Lambda] = Heat_Transfer_Procedure_Yarn(Work_Directory, CPU, NODES, ELEMENTS, Material_Properties, Coordinate_Systems, RVE_Coordinates, Node_Merge_Tolerance, V_resin);
        
        if (ErrCode < 0)
            return;
        end
        
        Yarn_Properties.Density = Density;
        
        Yarn_Properties.Cp = Cp;
        
        Yarn_Properties.Lambda_11 = Lambda(1,1);
        Yarn_Properties.Lambda_22 = Lambda(2,2);
        Yarn_Properties.Lambda_33 = Lambda(3,3);
    
end
